covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_adaptation_sens <- causalsens(climchg_adaptation_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Harris-Public-Policy Dropbox/Christopher Blair/Climate and Laws/Hurricanes/Figures/sens_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Adaptation Policies")
dev.off()
png(file="/Users/cb2257/Harris-Public-Policy Dropbox/Christopher Blair/Climate and Laws/Hurricanes/Figures/sensr_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Adaptation Policies")
dev.off()
#SCIENCE OF CLIMATE CHANGE
climchg_science_lm <- lm(climchg_science ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_science_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_science_sens <- causalsens(climchg_science_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Harris-Public-Policy Dropbox/Christopher Blair/Climate and Laws/Hurricanes/Figures/sens_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Science of Climate Change")
dev.off()
png(file="/Users/cb2257/Harris-Public-Policy Dropbox/Christopher Blair/Climate and Laws/Hurricanes/Figures/sensr_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Science of Climate Change")
dev.off()
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(maptools)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(maptools)
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(sp)
library(raster)
wbluescale <- c("#ffffff", "#08519c")
worangescale <- c("#ffffff", "#fc9272")
wgreenscale <- c("#ffffff", "#006d2c")
wgrayscale <- c("#ffffff", "#AAAAAA")
wpinkscale <- c("#ffffff", "#FF66B2")
wpurpscale <- c("#ffffff", "#9e9ac8")
wredscale <- c("#ffffff", "#de2d26")
bwhitescale6 <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#ffffff")
wbluescale6 <- c("#ffffff", "#c6dbef", "#9ecae1", "#6baed6", "#3182bd", "#08519c")
setwd("/Users/cb2257/Harris-Public-Policy Dropbox/Christopher Blair/Climate and Laws/Hurricanes/")
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
library(maptools)
library(sp)
library(raster)
wbluescale <- c("#ffffff", "#08519c")
worangescale <- c("#ffffff", "#fc9272")
wgreenscale <- c("#ffffff", "#006d2c")
wgrayscale <- c("#ffffff", "#AAAAAA")
wpinkscale <- c("#ffffff", "#FF66B2")
wpurpscale <- c("#ffffff", "#9e9ac8")
wredscale <- c("#ffffff", "#de2d26")
bwhitescale6 <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#ffffff")
wbluescale6 <- c("#ffffff", "#c6dbef", "#9ecae1", "#6baed6", "#3182bd", "#08519c")
setwd("/Users/cb2257/Desktop/APSR Replication/")
ian <- st_read("Shapefiles/hurricane_sample/hurricane_sample.shp")
rm(list=ls())
library(terra)
library(tmap)
library(tmaptools)
library(sf)
library(leaflet)
library(haven)
library(rgdal)
library(maptools)
library(sp)
library(raster)
wbluescale <- c("#ffffff", "#08519c")
worangescale <- c("#ffffff", "#fc9272")
wgreenscale <- c("#ffffff", "#006d2c")
wgrayscale <- c("#ffffff", "#AAAAAA")
wpinkscale <- c("#ffffff", "#FF66B2")
wpurpscale <- c("#ffffff", "#9e9ac8")
wredscale <- c("#ffffff", "#de2d26")
bwhitescale6 <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#ffffff")
wbluescale6 <- c("#ffffff", "#c6dbef", "#9ecae1", "#6baed6", "#3182bd", "#08519c")
setwd("/Users/cb2257/Desktop/APSR Replication/")
ian <- st_read("Shapefiles/hurricane_sample/hurricane_sample.shp")
st_is_valid(ian, reason = TRUE)
ian <- st_make_valid(ian)
south <- st_read("Shapefiles/southern_states/southern_states.shp")
st_is_valid(south, reason = TRUE)
south <- st_make_valid(south)
eyepath <- st_read("Shapefiles/ian_eyepath/ian_eyepath.shp")
st_is_valid(eyepath, reason = TRUE)
eyepath <- st_make_valid(eyepath)
ida <- st_read("Shapefiles/ida_eyepath/ida_eyepath.shp")
st_is_valid(ida, reason = TRUE)
ida <- st_make_valid(ida)
nicole <- st_read("Shapefiles/nicole_eyepath/nicole_eyepath.shp")
st_is_valid(nicole, reason = TRUE)
nicole <- st_make_valid(nicole)
county <- read_dta(file="Data/Survey/map_descriptives.dta")
joined <- merge(ian, county, by.x="GEOID", by.y="geoid")
# Map of Respondents -- Figure 1 in Main Paper
respondents <- tm_shape(joined)+
tm_borders()+
tm_fill("response", title = "# of Respondents", palette=wgreenscale, style="fixed",
breaks=c(0, 1, 4, 11, 27, 280), interval.closure="left",
labels=c("0", "1-3", "4-10", "11-26", "27+"))+
tm_layout(" ", title.size=1, frame = F)
respeye<-respondents+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+respeye
tmap_save(combined, "Figures/respondent_map.png")
# Map of Respondents2 -- Figure A-12 in Main Appendix
respondents2 <- tm_shape(joined)+
tm_borders()+
tm_fill("response2", title = "# of Respondents", palette=wgreenscale, style="fixed",
breaks=c(0, 1, 4, 11, 27, 58), interval.closure="left",
labels=c("0", "1-3", "4-10", "11-26", "27+"))+
tm_layout(" ", title.size=1, frame = F)
respeye<-respondents2+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+respeye
tmap_save(combined, "Figures/respondent2_map.png")
# Map of Windswath -- Figure 2c in Main Paper
windswath <- tm_shape(joined)+
tm_borders()+
tm_fill("ian_windswath", title = "Windswath", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4), interval.closure="left",
labels=c("Sub-Cyclonic Winds (0-33 Knots)", "Tropical Storm Winds (34+ Knots)", "Violent Storm Winds (50+ Knots)", "Hurricane Winds (64+ Knots)"))+
tm_layout(" ", title.size=1, frame = F)
windeye<-windswath+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+windeye
tmap_save(combined, "Figures/windswath_map.png")
# Map of FEMA Aid -- Figure 2f in Main Paper
fema <- tm_shape(joined)+
tm_borders()+
tm_fill("fema_aid", title = "FEMA Assistance", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4, 5), interval.closure="left",
labels=c("None", "Public (Category B)", "Public (Category A-G)", "Individual and Public (Category B)", "Individual and Public (Category A-G)"))+
tm_layout(" ", title.size=1, frame = F)
femaeye<-fema+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+femaeye
tmap_save(combined, "Figures/fema_map.png")
# Map of Distance -- Figure 2d in Main Paper
distance <- tm_shape(joined)+
tm_borders()+
tm_fill("ian_distance_mi", title = "Distance to Eyepath", palette=bwhitescale6, style="fixed",
breaks=c(0, 77.05172, 190.0623, 895.7505, 1154.872, 1393.435, 2000), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
distanceye<-distance+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+distanceye
tmap_save(combined, "Figures/distance_map.png")
# Map of Evacuations -- Not presented in paper
evacuation <- tm_shape(joined)+
tm_borders()+
tm_fill("ian_evacuation", title = "Evacuation Order", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2), interval.closure="left",
labels=c("No", "Yes"))+
tm_layout(" ", title.size=1, frame = F)
evaceye<-evacuation+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+evaceye
tmap_save(combined, "Figures/evacuation_map.png")
# Map of Deaths -- Not presented in paper
deaths <- tm_shape(joined)+
tm_borders()+
tm_fill("i_ian_deaths", title = "Fatalities", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2), interval.closure="left",
labels=c("No", "Yes"))+
tm_layout(" ", title.size=1, frame = F)
deathseye<-deaths+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+deathseye
tmap_save(combined, "Figures/i_deaths_map.png")
# Map of Deaths -- Figure 2e in Main Paper
deaths <- tm_shape(joined)+
tm_borders()+
tm_fill("ian_deaths", title = "# of Fatalities", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 6, 10, 73), interval.closure="left",
labels=c("0", "1", "2-5", "6-9", "10+"))+
tm_layout(" ", title.size=1, frame = F)
deathseye<-deaths+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+deathseye
tmap_save(combined, "Figures/deaths_map.png")
# Map of Storm Surge -- Figure 2b in Main Paper
surge <- tm_shape(joined)+
tm_borders()+
tm_fill("ian_surge", title = "Storm Surge", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4, 5), interval.closure="left",
labels=c("None", "1-3 Feet", "3-6 Feet", "6-9 Feet", "9+ Feet"))+
tm_layout(" ", title.size=1, frame = F)
surgeeye<-surge+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+surgeeye
tmap_save(combined, "Figures/surge_map.png")
# Map of Index -- Figure 2a in Main Paper
index <- tm_shape(joined)+
tm_borders()+
tm_fill("hurricane_pctile", title = "Hurricane Exposure Index", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4, 5, 6), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
indexeye<-index+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+indexeye
tmap_save(combined, "Figures/index_map.png")
# Map of Binary Index -- Figure A-5 in Main Appendix
binary <- tm_shape(joined)+
tm_borders()+
tm_fill("i_index", title = "Hurricane Exposure", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2), interval.closure="left",
labels=c("No", "Yes"))+
tm_layout(" ", title.size=1, frame = F)
binaryeye<-binary+tm_shape(eyepath)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+binaryeye
tmap_save(combined, "Figures/i_index_map.png")
# Map of Ida -- Figure A-9 in Main Appendix
idaindex <- tm_shape(joined)+
tm_borders()+
tm_fill("ida_pctile", title = "Hurricane Exposure Index", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4, 5, 6), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
idaeye<-idaindex+tm_shape(ida)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+idaeye
tmap_save(combined, "Figures/ida_map.png")
# Map of Nicole  -- Figure SI-3 in Additional Supplement
nicoleindex <- tm_shape(joined)+
tm_borders()+
tm_fill("nicole_pctile", title = "Hurricane Exposure Index", palette=wbluescale, style="fixed",
breaks=c(0, 1, 2, 3, 4, 5, 6), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
nicoleeye<-nicoleindex+tm_shape(nicole)+tm_lines(col = "red", lwd = 2, lty = "dashed")
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+nicoleeye
tmap_save(combined, "Figures/nicole_map.png")
# Map of 2020 Republican Voteshare -- Not presented in paper
repvote <- tm_shape(joined)+
tm_borders()+
tm_fill("repvote20", title = "2020 Republican Voteshare", palette=wredscale, style="fixed",
breaks=c(0, .2, .4, .6, .8, 1), interval.closure="left",
labels=c("< 20%", "20-39%", "40-59%", "60-79%", "> 80%"))+
tm_layout(" ", title.size=1, frame = F)
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+repvote
tmap_save(combined, "Figures/repvote_map.png")
repwon <- tm_shape(joined)+
tm_borders()+
tm_fill("repvote20", title = "2020 Presidential Election", palette=wredscale, style="fixed",
breaks=c(0, .5, 1), interval.closure="left",
labels=c("Biden Won", "Trump Won"))+
tm_layout(" ", title.size=1, frame = F)
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+repwon
tmap_save(combined, "Figures/repwon_map.png")
# Map of Relief Aid -- Not presented in paper
aid <- tm_shape(joined)+
tm_borders()+
tm_fill("ihpAmount", title = "FEMA Individual \n Assistance per Applicant", palette=wpinkscale, style="fixed",
breaks=c(0, .00000001, 495, 766, 1700, 2453, 2958), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+aid
tmap_save(combined, "Figures/femaind_map.png")
aid <- tm_shape(joined)+
tm_borders()+
tm_fill("numberobligatedprojects", title = "FEMA Public \n Assistance Projects", palette=wpinkscale, style="fixed",
breaks=c(0, 1, 2, 8, 12, 21, 84), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+aid
tmap_save(combined, "Figures/femaproj_map.png")
aid <- tm_shape(joined)+
tm_borders()+
tm_fill("obligated", title = "FEMA Public \n Assistance per Project", palette=wpinkscale, style="fixed",
breaks=c(0, .00000001, 15077, 47303, 150685, 639896, 3729566), interval.closure="left",
labels=c("< 10th Percentile", "10-25th Percentile", "25-50th Percentile", "50-75th Percentile", "75-90th Percentile", "> 90th Percentile"))+
tm_layout(" ", title.size=1, frame = F)
combined<-tm_shape(south)+tm_fill(col = "gray83")+tm_borders(col = "black", lwd = 1, lty = "solid")+aid
tmap_save(combined, "Figures/femapub_map.png")
remove(list=ls())
setwd("/Users/cb2257/Desktop/APSR Replication/")
library("lme4")
library("bayestestR")
library("tidyverse")
library("dplyr")
library("leaflet")
library("haven")
library("leaflet")
library("ggplot2")
library("mosaic")
library("estimatr")
data <- read_dta("Data/Survey/survey_final.dta")
estsample <- data[which(data$estsample==1),]
estsample <- estsample %>% mutate(std_age = scale(age))
estsample <- estsample %>% mutate(std_ideology = scale(ideology))
varlist <- c("republican", "democrat", "woman", "highschool", "college", "std_age", "std_ideology", "race_white", "race_black", "race_latinx", "employed", "native", "religiosity", "empathy", "pol_interest")
#Perform Equivalence Tests -- Figure A-2 in Main Appendix
est_equiv <- function(varname) {
require(parameters)
require(stats)
#Define Empirical Model
f <- paste0("hurricane_index1 ~ factor(start_date) + factor(statefp) + ", varname)
#Estimate Model With Weights
m.wt <- lm(f, data = estsample, weights = eweight)
#Conduct Equivalence Test
equiv.wt <- equivalence_test(m.wt, rule = "classic", ci=.975, range(c(-.15,.15)))
##Extract Model Information
df.wt <- data.frame(equiv.wt)
df.wt$sample <- " "
df.wt$estimate <- coef(m.wt)
#subset intercept since not of interest
df.out <- subset(df.wt, Parameter!=c("(Intercept)", "factor(start_date)2022-08-22" , "factor(start_date)2022-08-23", "factor(start_date)2022-08-24", "factor(start_date)2022-08-26", "factor(start_date)2022-08-27", "factor(start_date)2022-08-29", "factor(start_date)2022-09-06", "factor(start_date)2022-09-19", "factor(start_date)2022-09-26", "factor(start_date)2022-09-27", "factor(start_date)2022-09-28", "factor(start_date)2022-09-29", "factor(start_date)2022-10-03", "factor(start_date)2022-10-04", "factor(start_date)2022-10-07", "factor(start_date)2022-10-12", "factor(start_date)2022-10-21", "factor(start_date)2022-10-27"))
return(df.out)
}
#loop over variable list
equiv.out <- lapply(varlist, est_equiv)
equiv.out <- do.call("rbind", equiv.out)
equiv.out
republican <- subset(equiv.out,Parameter=="republican")
democrat <- subset(equiv.out,Parameter=="democrat")
highschool <- subset(equiv.out,Parameter=="highschool")
college <- subset(equiv.out,Parameter=="college")
woman <- subset(equiv.out,Parameter=="woman")
age <- subset(equiv.out,Parameter=="std_age")
ideology <- subset(equiv.out,Parameter=="std_ideology")
white <- subset(equiv.out,Parameter=="race_white")
black <- subset(equiv.out,Parameter=="race_black")
latinx <- subset(equiv.out,Parameter=="race_latinx")
employed <- subset(equiv.out,Parameter=="employed")
native <- subset(equiv.out,Parameter=="native")
religiosity <- subset(equiv.out,Parameter=="religiosity")
empathy <- subset(equiv.out,Parameter=="empathy")
interest <- subset(equiv.out,Parameter=="pol_interest")
equiv.out<-rbind(republican, democrat, highschool, college, woman, age, ideology, white, black, latinx, employed, native, religiosity, empathy, interest)
equiv.out$Parameter  = factor(equiv.out$Parameter, levels=c("pol_interest", "empathy", "religiosity", "native", "employed", "race_latinx", "race_black", "race_white", "std_ideology", "std_age", "college", "highschool", "woman", "democrat", "republican"))
equiv.out$Estimate <- equiv.out$sample
equiv.out$Hurricane <- equiv.out$estimate
renames <- c("republican" = "Republican", "democrat" = "Democrat", "highschool" = "High School Graduate", "college" = "College Graduate", "woman" = "Woman", "std_age" = "Age", "std_ideology" = "Political Ideology", "race_white" = "White", "race_black" = "Black", "race_latinx" = "Latinx", "employed" = "Employed", "native" = "Native Born", "religiosity" = "Religiosity", "empathy" = "Empathy", "pol_interest" = "Political Interest")
equiv.out <- equiv.out %>% mutate(Parameter = str_replace_all(Parameter, renames))
#create balance plot
rope_low <- equiv.out$ROPE_low[[1]]
rope_high <- equiv.out$ROPE_high[[1]]
equiv.out %>%
ggplot() +
geom_hline(yintercept=rope_low, lty="dashed") +
geom_hline(yintercept=rope_high, lty="dashed") +
geom_pointrange(aes(x=Parameter, y=Hurricane, ymin=CI_low, ymax=CI_high, color=Estimate),
position=position_dodge(.5)) +
coord_flip() +
scale_y_continuous(name="Hurricane Exposure", limits = c(-.3, .3), breaks=c(-.3,-.2,-.1,0,.1,.2,.3))+
theme(legend.position="none")+
theme_bw(base_size = 14)
ggsave(file="equivalency.png", path="Figures/", device='png', dpi=300)
rm(list = ls())
library(causalsens)
library(DirectEffects)
library(haven)
library(foreign)
library(rio)
library(tidyverse)
library(sjlabelled)
library(memisc)
library(expss)
library(dplyr)
library(tidyverse)
library(magrittr)
library(sandwich)
library(miceadds)
library(stargazer)
ian <- import("/Users/cb2257/Desktop/APSR Replication/Data/Survey/survey_final.dta")
estsample <- subset(ian, estsample == 1)
core <- c("geoid", "start_date")
estsample$treatment <- estsample$hurricane_index1 * estsample$post
estsample[core] <- lapply(estsample[core], factor)
#BLACKWELL (2014) SENSITIVITY TESTS
#CLIMATE MIGRATION ISSUE IMPORTANCE -- Figure A-8 in Main Appendix
climmig_importance_lm <- lm(climmig_importance ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climmig_importance_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
pmodel <- glm(treatment ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_pmodel <- summary(pmodel, cluster = "geoid")
summary_pmodel
alpha <- seq(-.15, .15, by = 0.01)
climmig_importance_sens <- causalsens(climmig_importance_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climmig_importance.png", width=6, height=4, units="in", res=1200)
plot(climmig_importance_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Migration Issue Importance")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climmig_importance.png", width=6, height=4, units="in", res=1200)
plot(climmig_importance_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Migration Issue Importance")
dev.off()
#CLIMATE MIGRATION POLICY ACTION -- Figure A-8 in Main Appendix
climmig_action_lm <- lm(climmig_action ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climmig_action_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climmig_action_sens <- causalsens(climmig_action_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climmig_action.png", width=6, height=4, units="in", res=1200)
plot(climmig_action_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Migration Policy Action")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climmig_action.png", width=6, height=4, units="in", res=1200)
plot(climmig_action_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Migration Policy Action")
dev.off()
#CLIMATE CHANGE ISSUE IMPORTANCE -- Figure A-8 in Main Appendix
climchg_importance_lm <- lm(climchg_importance ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_importance_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_importance_sens <- causalsens(climchg_importance_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_importance.png", width=6, height=4, units="in", res=1200)
plot(climchg_importance_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Issue Importance")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_importance.png", width=6, height=4, units="in", res=1200)
plot(climchg_importance_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Issue Importance")
dev.off()
#CLIMATE CHANGE POLICY ACTION -- Figure A-8 in Main Appendix
climchg_action_lm <- lm(climchg_action ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_action_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_action_sens <- causalsens(climchg_action_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_action.png", width=6, height=4, units="in", res=1200)
plot(climchg_action_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Policy Action")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_action.png", width=6, height=4, units="in", res=1200)
plot(climchg_action_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Policy Action")
dev.off()
#CLIMATE CHANGE MITIGATION POLICIES -- Figure A-8 in Main Appendix
climchg_mitigation_lm <- lm(climchg_mitigation ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_mitigation_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_mitigation_sens <- causalsens(climchg_mitigation_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_mitigation.png", width=6, height=4, units="in", res=1200)
plot(climchg_mitigation_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Mitigation Policies")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_mitigation.png", width=6, height=4, units="in", res=1200)
plot(climchg_mitigation_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Mitigation Policies")
dev.off()
#CLIMATE CHANGE ADAPTATION POLICIES -- Figure A-8 in Main Appendix
climchg_adaptation_lm <- lm(climchg_adaptation ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_adaptation_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_adaptation_sens <- causalsens(climchg_adaptation_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Adaptation Policies")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Adaptation Policies")
dev.off()
#SCIENCE OF CLIMATE CHANGE -- Figure A-8 in Main Appendix
climchg_science_lm <- lm(climchg_science ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_science_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary
alpha <- seq(-.15, .15, by = 0.01)
climchg_science_sens <- causalsens(climchg_science_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Science of Climate Change")
dev.off()
png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Science of Climate Change")
dev.off()
